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1  INTRODUCTION 


SHADO  is  a  new  module  for  use  in  computing  the  particle  wake  behind 
an  object  in  Low  Earth  Orbit  (LEO)  traveling  at  mesosonic  speeds.  At 
mesosonic  velocities,  a  spacecraft  will  sweep  away  ambient  ions  since  it  is 
moving  faster  that  the  ion  thermal  velocity  and  much  slower  than  the  elec¬ 
tron  thermal  velocity.  This  condition  is  encountered  in  LEO  where  the  night 
side  plasma  can  have  a  Debye  length  of  2  cm  while  the  ion  Debye  length  is 
roughly  one  fifth  the  spacecraft  dimension.  The  motion  of  ions  to  fill  in  the 
electron  rich  wake  behind  the  spacecraft  is  quickly  dominated  by  the  electric 
field  in  the  wake  which  is  obtained  by  solving  Poisson’s  equation  [1], 

Two  simplifying  assumptions  are  made  in  computing  the  ion  densities  in 
the  wake.  The  first  assumption  is  the  neutral  ion  approximation  which  as¬ 
sumes  ions  travel  in  straight  lines  and  are  stopped  if  they  impact  the  space¬ 
craft.  This  neglects  electric  field  effects  on  ion  trajectories.  The  second 
assumption  is  that  ion  filling  of  the  wake  is  due  to  electric  field  accelera¬ 
tion,  thereby  neglecting  ion  momentum.  The  neutral  ion  approximation  is 
expressed  by  the  following  equation: 

/.(*,»)  =  g(x,n)  fioi V)  (1) 

where  /,(x,  v)  is  the  ion  distribution  function  at  a  point  x  in  space  for  a 
velocity  v  and  fto{v)  is  the  unperturbed  velocity  distribution  function  for 
a  drifting  Maxwellian.  The  function  g{x,fl)  has  a  value  of  zero  if  a  ray 
starting  from  x  in  direction  fl  would  strike  the  spacecraft;  otherwise  it  is 
one.  The  charge  density  is  obtained  by  integration  over  velocities: 


".(*)  =  f  f(x,v)  =  fx0{v^)v2di/ 
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(2) 


The  object  being  studied  is  described  by  a  collection  of  surface  elements 
which  are  referred  to  as  object  definition  surface  elements.  In  order  to  deter¬ 
mine  the  charge  density  at  a  point  in  space,  it  is  necessary  to  determine  the 
solid  angle  which  is  occupied  by  the  object  and  performing  the  integration 
in  equation  2.  The  problem  is  reduced  to  finding  g( x,fl)  given  x  due  to  the 
neutral  ion  approximation. 

The  older  algorithm  computed  the  charge  density  for  a  given  point  in 
space  by  determining  the  solid  angle  obstructed  by  each  surface  in  the  object 


definition  obstructed  when  viewed  from  x.  This  was  accomplished  by  defin¬ 
ing  a  spherical  grid  centered  at  f  whose  9  =  0,0  =  0  ray  pointed  directly 
into  the  oncoming  ions.  The  grid  was  subdivided  into  1  degree  elements 
in  the  polar  angle  6,  and  10  degree  elements  in  <t>.  The  greater  resolution 
in  9  was  required  since  the  function  /„  changes  rapidly  as  the  polar  angle 
increases.  Each  surface  was  then  projected  onto  this  spherical  grid  and  a 
number  of  extra  vertices  added  to  account  for  the  curvature  induced  by  the 
projection  of  straight  lines  onto  a  curved  surface.  Each  element  in  the  sphere 
which  was  covered  by  the  projected  surface  element  was  then  given  a  value 
of  9s,*  —  0  and  integrated  over  9  and  <j>  to  obtain  the  charge  density.  [2] 

While  this  approach  was  practical  for  simple  objects,  it  became  clear 
that  complex  objects  with  many  surface  elements  would  be  increasingly 
difficult  to  handle.  This  is  especially  true  when  several  million  x  coordinates 
are  required  to  adequately  compute  the  particle  wake  behind  a  spacecraft. 
Another  application  for  which  this  approach  is  inadequate  is  problem  of 
sheath  shadowing.  The  sheath  around  the  moving  object  can  be  described 
using  surface  elements  and  treated  as  an  object.  However,  very  many  surface 
elements  are  required  and  the  old  algorithm  is  too  cumbersome  to  compute 
both  the  particle  wake  behind  the  object  and  the  sheath  shadowing. 

A  new  approach  was  then  developed  which  greatly  simplifies  the  object 
definition  so  that  the  charge  density  may  be  computed  much  more  rapidly. 
The  new  algorithm  also  relies  heavily  on  the  fact  that  the  ion  distribution, 
relative  to  the  spacecraft,  is  heavily  weighted  towards  the  ram  direction  at 
mesosonic  velocities.  Since  is  azimuthally  symmetric,  equation  2  can  be 
reduced  to: 

G/(*)  =  Y^fcum{9,)  6e>  (3) 

where  the  G/’s  are  called  ‘geometrical  ions’  obtained  by  dividing  equation  2 
by  the  unperturbed  ion  density  (n0, ).  If  a  point  x  is  completly  blocked  by  the 
spacecraft,  GI  will  be  zero,  while  GI  will  be  unity  far  from  the  spacecraft. 

The  function  fcum(9)  is  obtained  by  integrating  f>0{9,  ©,  |  v] )  over  <p  since 
it  is  azimuthally  symmetric,  and  normalizing  it  as  follows. 

fcum(9)  =  22fto(0„<t>„  v)Ac,  (4) 

£/cum(0)A0,  =  1  (5) 

s. 
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fcum  is  calculated  in  advance,  given  u,  and  is  stored  internally  as  an  array 
which  rapidly  decreases  as  9  increases  for  mesosonic  velocties.  The  value  of 
9  at  which  fcum  drops  to  0.01%  of  its  value  in  the  ram  direction  {9  =  0) 
is  computed  and  saved  as  #99%.  As  |u]  increases.  099%  decreases  indicating 
that  the  Maxwellian  distribution  of  ions  is  shifted  towards  the  ram  direction. 
Surfaces  lying  within  0  <  9  <  099%  from  the  perspective  of  x  will  therefore 
have  a  significant  effect  on  the  value  of  GI(x).  while  surfaces  positioned 
beyond  #99%  will  have  a  negligible  impact. 

Simplification  of  the  object  definition  is  achieved  by  first  projecting  the 
object  onto  a  surface  which  is  normal  to  the  ram  direction  and  then  overlay¬ 
ing  this  projection  with  a  new  two  dimensional  grid  as  illustrated  in  Figure 
4.  The  nodes  on  this  grid  are  arranged  in  an  equilaterallv  spaced  staggered 
array  resulting  in  hexagonal  elements.  These  hexagonal  elements  provide 
greater  angular  resolution  of  the  object  projection  than  an  equivalent  carte¬ 
sian  grid  system.  Each  surface  element  in  the  object  definition  is  then  se¬ 
quentially  projected  onto  this  hexagonal  grid  and  a  depth  is  computed  for 
each  grid  node  which  falls  inside  the  surface  element  projection.  This  depth 
represents  the  distance  from  the  projection  surface  to  the  surface  element 
along  a  ray  parallel  to  the  ram  direction.  The  projection  surface  is  placed  so 
that  the  forward-most  point  on  the  object  corresponds  to  zero  depth.  All  of 
these  calculations  are  performed  only  once  prior  to  evaluation  of  any  GTs. 

The  number  of  nodes  making  up  the  hexagonal  mesh  is  controlled  by  a 
parameter  .Y/,ex.  The  nodes  are  distributed  uniformly  over  the  rectangular 
area  on  the  projection  plane  which  just  encompasses  the  entire  object  pro¬ 
jection.  As  each  object  definition  surface  is  processed,  a  table  of  depths  for 
each  node  in  the  hexagonal  mesh  is  maintained.  After  all  the  surfaces  are 
projected,  the  table  for  each  node  is  sorted  by  increasing  depth.  Thus,  at 
each  node  on  the  hexagonal  grid  there  is  a  table  of  depths  representing  the 
intersection  of  a  ray  beginning  at  the  node  and  any  of  the  object  definition 
surfaces  in  the  path  along  the  anti-ram  direction.  Given  a  x  coordinate,  it  is 
then  simple  to  compute  its  distance  from  the  projection  plane  and  its  loca¬ 
tion  on  the  hexagonal  grid.  Comparing  the  depth  of  x  to  the  table  of  depths 
of  the  adjacent  hexagonal  grid  nodes  it  is  quickly  determined  whether  x  is 
in  front  of.  behind,  inside,  or  adjacent  to  the  object. 

If  the  reference  point,  f.  is  in  front  of  the  object,  then  GI(x)  will  he 
unity  because  all  angles  9  pointing  to  the  object  definition  surfaces  will  be 
much  larger  than  999%.  If  the  reference  point  is  sufficiently  far  away  from 


the  object  such  that  the  entire  object  lies  outside  of  the  9 99%  limit,  then 
GJ  is  set  to  unity  as  well.  If  the  reference  point  is  inside  the  object,  then 
G/  is  set  to  zero.  Since  many  of  the  points  in  the  object  mesh  fall  into  one 
of  these  categories,  little  computational  effort  is  expended  in  determining 
the  GJ(x)’s.  For  the  remaining  points  in  the  object  mesh,  it  is  necessary  to 
determine  where  the  object  is  relative  to  x  and  determine  how  much  of  the 
view  towards  the  am  direction  is  blocked. 

Much  of  the  computational  speed  improvement  in  the  SHADO  module 
is  attributable  to  the  simplification  of  the  object  representation  already  de¬ 
scribed.  Given  an  x  coordinate,  say  directly  behind  the  object,  a  scan  along 
a  number  of  rays  controlled  by  a  parameter  iVp^,,  is  initiated.  Each  ray 
begins  at  x  and  proceeds  along  a  fixed  angle  in  the  plane  of  the  projection 
surface.  At  intervals  just  less  than  the  hexagonal  mesh  element  size,  the  sur¬ 
rounding  hexagonal  nodes  are  tested  to  determine  whether  they  are  over  the 
object  or  not.  If  the  table  of  depths  for  one  or  more  of  these  nodes  is  empty, 
then  this  indicates  that  the  scan  has  proceeded  beyond  the  projection  of 
the  object.  A  quick  calculation  is  made  which  approximates  the  distance  to 
the  object  outline  (accurate  to  half  the  hexagonal  mesh  dimension)  and  the 
depth  to  the  object  at  this  point  is  computed.  This  gives  the  angle  9  from 
x  to  the  boundary  of  the  object  and  the  contribution  to  G/(x)  along  this 
value  of  <f>  is  determined.  Summing  for  each  of  the  Nphi  values  of  <*>  gives 
the  value  the  geometrical  ion  at  x. 

A  number  of  variations  on  the  scheme  outlined  above  are  possible  based 
on  whether  the  point  x  is  directly  behind  (the  table  of  depths  on  all  nodes 
adjacent  to  x  are  non-zero),  or  off  to  one  side  of  the  object.  The  same  basic 
procedure  is  followed  however,  by  scanning  in  rays  along  the  projection 
plane,  determining  where  the  object  is  relative  to  x  on  this  plane,  and  then 
determining  the  angles,  9,  locating  the  boundaries  of  the  object  in  three 
dimensions  relative  to  x.  The  accuracy  of  the  object  description  using  this 
algorithm  is  controlled  by  the  number  of  hexagonal  grid  nodes  (A \ez)  placed 
on  the  hexagonal  grid  and  the  number  of  scan  angles  ( A’p/l<)  used  when 
locating  the  object  relative  to  x. 

There  are  also  a  number  of  shortcuts  used  in  the  implementation  of  this 
algorithm.  The  major  one  being  that  the  hexagonal  grid  is  skewed  such 
that  the  y  axis  is  at  a  60  degree  angle  with  the  x  axis  on  the  projection 
plane.  This  results  in  a  cartesian  spacing  of  the  grid  nodes,  which  when 
normalized  to  the  distance  between  grid  points,  puts  the  skewed  grid  nodes 
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on  integer  coordinates.  During  the  SHADO  setup,  the  object  definition 
surface  elements  are  actually  projected  onto  this  skewed  coordinate  system, 
and  all  of  the  scans  along  lines  of  constant  <f>  take  place  in  this  coordinate 
system  as  well.  In  the  example  problems  which  follow,  the  objects  are  shown 
projected  on  this  skewed  mesh  with  a  resulting  distortion  of  the  object. 

2  SAMPLE  CALCULATIONS 

A  number  of  test  cases  have  been  run  to  determine  the  accuracy  limits, 
the  computational  speed  improvements,  and  to  compare  results  against  the 
older  algorithm.  In  particular,  comparisons  were  made  for  a  24  surface 
quasisphere,  a  200  surface  shuttle  orbiter,  and  a  1700  surface  representation 
of  the  National  Aerospaceplane.  The  quasisphere  was  used  to  investigate 
accuracy  limits  for  varying  values  of  N^tx  (the  hexagonal  grid  resolution) 
and  JV*  (the  angular  resolution  in  the  hexagonal  grid).  The  shuttle  orbiter 
model  was  to  demonstrate  the  speed  improvements  over  the  older  algorithm, 
and  the  National  Aerospaceplane  model  served  as  a  demonstration  of  the 
flexibility  of  this  new  algorithm  for  problems  of  great  complexity. 

-  2.1  Quasisphere  Results 

The  quasisphere  is  a  simple  object  made  up  of  24  surfaces  which  approximate 
a  sphere  with  a  diameter  of  12.  Calculations  for  both  the  new  and  old 
algorithms  consisted  of  calculating  the  density  along  a  line  5  radii  behind 
the  sphere  and  perpendicular  to  the  ram  direction.  The  line  begins  roughly 
one  diameter  off  the  axis  of  symmetry  and  extends  to  the  axis  of  symmetry. 
Figure  5  shows  the  results  obtained  for  this  problem  from  the  old  algorithm. 
The  curve  represents  the  drop  in  GI  as  a  traverse  is  made  from  20  units 
off  the  centerline  to  the  centerline.  The  quasisphere  is  centered  along  the 
x=10  value  on  the  figure.  The  Mach  number  used  was  8.1  and  showed  good 
agreement  with  the  experimental  results  of  Fournier  and  Pigache  [l]. 

Figure  6  shows  the  projection  of  the  quasisphere  on  the  hexagonal  mesh 
used  in  the  new  algorithm.  The  grid  size  is  a  function  of  the  Nhei  parameter 
and  is  set  to  1000  which  means  that  the  hex  grid  is  scaled  such  that  no  more 
than  1000  elements  are  used  to  contain  the  projection  of  the  object.  De¬ 
pending  on  the  aspect  ratio  of  the  object  projection,  the  number  of  elements 
will  generally  be  less  than  N^ex  by  a  small  amount.  The  stars  in  figure  6 


represent  the  vertices  of  a  three  dimensional  box  which  just  contains  the 
object.  Due  to  the  direct  alignment  of  the  quasisphere  with  respect  to  the 
ram  direction,  four  of  the  eight  vertices  of  this  box  are  stacked  directly  on 
top  of  the  other  four. 

Figures  7-10,  show  the  dependency  of  the  solution  on  the  angular  reso¬ 
lution  with  N#  set  to  8,  16,  32  and  64  respectively.  The  effect  of  increasing 
angular  resolution  is  a  smoothing  out  of  the  density  contour.  The  bump  in 
figure  7  and  to  a  lesser  extent  in  figure  8,  is  caused  by  the  interception  of  a 
single  scan  angle  when  it  first  finds  that  it  is  blocked  by  the  object. 

Another  study  of  the  sensitivity  to  the  A\ex  parameter  for  grid  resolution 
showed  no  improvements  for  values  higher  than  the  value  of  1000  shown  in 
figure  6,  although  the  density  contour  roughness  was  slightly  reduced  as  A \ez 
was  increased.  A  general  rule  of  thumb  would  be  to  set  .\heT  to  adequately 
resolve  the  smallest  feature  of  interest  on  the  object. 

2.2  Shuttle  Orbiter  Results 

Id  order  to  compare  the  speed  of  execution  for  both  the  new  and  the  old 
algorithms,  a  3-D  calculation  of  the  ion  density  around  the  shuttle  orbiter 
was  run.  Figure  11  shows  the  shuttle  object  definition  and  its  projection  on 
the  hexagonal  mesh.  In  the  wake  calculation  shown  in  Figure  12.  the  orbiter 
is  oriented  with  its  bay  doors  open  and  in  the  ram  direction.  The  value  of 
A^ei  was  1000  and  the  value  for  was  16.  Figure  8  shows  the  resulting 
particle  wake  through  a  section  taken  from  approximately  midships  (the  nose 
of  the  orbiter  points  into  the  page,  the  wings  are  vertical).  The  coarseness  in 
the  density  contours  is  due  to  the  low  resolution  used  in  selecting  f  values. 
The  Mach  number  was  8.0  for  this  problem  and  A'0  was  set  to  16. 

The  speed  of  the  SHADO  algorithm  was  clearly  demonstrated  in  this 
example.  The  old  algorithm  required  roughly  30  times  more  CPU  time  to 
calculate  ion  density  for  the  roughly  100,000  points  than  did  the  SHADO 
algorithm.  On  a  Ridge  32  computer,  the  computation  time  was  reduced 
from  approximately  30  hours  to  under  one  hour. 


2.3  National  Aerospaceplane  Results 

As  a  demonstration  of  the  versatility  of  this  new  algorithm,  a  particle  wake 
was  computed  for  a  detailed  rendition  of  the  National  Aerospaceplane.  de¬ 
scribed  by  over  1700  surface  elements.  The  angle  of  attack  was  20  degrees 
at  Mach  8.  Figure  13  shows  the  object  definition  used  and  the  projection 
onto  the  hexagonal  mesh.  Figure  14  shows  the  particle  wake  behind  the 
aerospaceplane  taken  though  the  midplane  and  includes  a  silhouette  of  the 
object.  This  calculation  required  A\er  to  be  set  at  5000  in  order  to  ade¬ 
quately  capture  the  particle  wake  at  the  sharp  nose.  A'0  was  set  to  64  and 
a  dense  spatial  resolution  in  x  was  used  in  order  to  obtain  smooth  density 
contours.  The  calculation  shown  in  the  figure  represents  over  10,000  spatial 
coordinates  and  required  a  little  over  7  minutes  of  CPU  time  on  a  Ridge  32 
computer. 

This  example  demonstrates  the  strength  of  this  new  algorithm  in  its 
ability  to  deal  with  complex  objects  in  an  efficient  manner.  The  SHDSET 
setup  routine  first  reduces  the  problem  from  a  collection  of  many  surfaces 
to  a  projected  surface  of  order  A\eI  elements.  Once  this  preprocessing  is 
accomplished,  the  algorithm  is  constrained  only  by  the  mesh  size  as  dictated 
by  the  A\et  parameter  and  the  angular  resolution  as  dictated  by  the  A'* 
parameter.  The  computational  efficiency  of  the  old  algorithm  at  each  point 
in  space  was  directly  linked  to  the  number  of  surfaces  used  in  the  object 
representation  which  often  required  that  complex  objects  be  simplified  by 
combining  surface  elements  into  larger  elements. 

3  SHADO  STRUCTURE 

The  SHADO  routines  are  organized  in  a  highly  modular  fashion.  The  fol 
lowing  sections  detail  the  organization  of  the  routines  and  describes  their 
function.  The  overall  structure  diagram  for  SHADO  is  shown  in  Figure  1. 
The  driver  is  the  POLAR  code.  SHDSET  is  the  routine  which  processes 
the  object  definition  surface  elements,  and  SHADO  is  the  routine  which 
computes  the  GI  given  x. 
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Figure  1:  Shado  Structure  Diagram 

3.1  SHDSET  -  Set  Up  Routine 

Initialization  for  the  SHADO  routine  takes  place  in  SHDSET  and  consists 
of  looping  over  each  surface  defining  an  object,  projecting  that  surface  onto 
the  plane  perpendicular  to  the  ram  direction,  and  determining  which  nodes 
on  a  hexagonal  grid  are  shaded  by  each  surface.  For  each  surface,  limits  are 
determined  which  define  a  box  in  the  hex  mesh  which  encloses  the  projected 
surface.  Only  points  inside  this  box  are  checked  for  shadowing  by  that 
surface. 

At  each  mesh  point  in  the  hexagonal  grid,  a  record  is  kept  of  the  height 
of  intersection  from  the  projection  surface  and  the  number  of  surfaces  which 
shade  a  given  mesh  point.  An  array  stores  the  z  values  of  each  surface 
intersection  sorted  by  ascending  node  number  and  then  by  ascending  z  value. 
Given  any  node,  it  is  then  simple  to  determine  if  this  node  is  shaded  in  the 
ram  direction  and  to  find  the  entry  and  exit  z  values. 

The  structure  diagram  for  SHDSET  is  shown  in  Figure  2. 
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Figure  2:  SHDSET  structure  diagram 

3.1.1  ROTATE 

Computes  the  rotation  matrix  required  to  convert  any  (x.y.z)  to  (x'.y'.z')  in 
which  the  x’-y’  plane  is  normal  to  the  Mach  vector.  The  x"  axis  is  taken  to 
be  in  the  x-z  plane. 

3.1.2  SPACEP  &  A2LOAD 

Reads  in  surface  element  data  for  the  National  aerospace  plane  using  the 
MSIO  library.  Alternatively,  A2  surfaces  may  be  read  using  routines  A2LOAD 
and  A2PREP.  These  routines  simply  read  in  data  from  MSIO  data  files  and 
load  an  array  of  vertices  for  each  surface.  These  routines  will  also  find  the 
minimum  and  maximum  x,y,  and  z  values  within  the  object  to  be  used  in 
defining  the  hexagonal  mesh.  Data  selection  is  determined  by  the  Ll’Sl'RF 
parameter  which  corresponds  to  the  logical  unit  number  for  the  surface  data 
19  for  the  A2’s,  9  for  the  spaceplane  data.  Surfaces  are  defined  by  coordi 
nates  for  surface  vertices  and  an  outward  pointing  normal  vector  Thi>  form 
is  exactly  how  A2  surfaces  are  defined,  for  the  spaceplane.  however,  the  mr 
mal  is  computed  using  three  distinct  vertices  to  compute  the  normal  vector 
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assuming  clockwise  node  numbering 


11.3  MSHDEF 

l ting  the  minimum  And  aiAumum  x.y  and  i  vah«  a  t>>>«  .»  define*}  »,tr 
which  the  object  is  fullj  contained  The  eight  vertices  mating  up  th.s  t. 
are  then  projected  onto  the  x'  >  plane  and  the  minimum  ana  max. m  in. 
and  v  values  are  found  A  hexagonal  mesh  ■»  then  established  »}.,<  i  t 
its  origin  at  iiwn  -  A  y.*,*  -  A  and  extends  up  to  -  A  i;  * 

where  A  is  the  mesh  spacing  A  is  « hosen  using  iteratioi.  -m  r,  tha' 
number  of  mesh  points  is  c  ,V4#I  The  mesh  is  then  sSf»«]  r*o  jegi.-*- 
give  a  square  mesh  and  normalired  b\  A  such  that  the  noo>>>  <n  th*  :t. 
are  integer  values 

3.14  FCAL 

Calculates  the  ion  densitv  factor  arras  for  a  point  .n  sp*.  »•  «Li  f.  .> 
stmeted  in  ail  directions  The  ion  deasits  factor  .»  computed  **  *  tun.  t. 
of  the  out  of  plane  angle  #  as  follows 


3.1.5  GF.TSRF 


Heads  in  the  surlaies  of  the  object  and  returns  the  internal  Limits  for  the 
hexagonal  mesh  points  which  must  be  checked  for  shadowing  Any  surface 
which  is  perpendicular  to  the  projection  surface  is  ignored  A  loop  inside 
SHDSFI  reads  m  the  number  of  surfaces  from  this  file,  and  GETSRF  gets 
a  sai.d  surta.  e  1  he  return  arguments  are  th»  coordinates  of  up  to  four 
o*n.i  e*  max.r.g  ip  *  >uria«  e  element  wtiich  is  not  perpendicular  to  the  x  .y 

piane 

3  1  «  INSIDE 

F  ills  an  arra>  whn  h  retains  all  of  the  surface  intersections  at  each  hexagonal 
grut  point  wh,<  h  is  shaded  b>  asurfaie  For  each  mesh  point  which  is  blocked 
Pc  '  rie  gi  .e:.  >urfai  th**  height  of  the  surface  above  the  projet  tion  surlai  e  is 
>  ompute.;  a:.d  *.t.>r»M  in  an  array  of  i  s  and  the  corresponding  node  number 
is  stored  ,n  a  lompanion  array  After  ail  of  the  surfaces  have  been  processed, 
these  tw«.  arrav s  are  sorted  bv  node  in  ascending  order  of  /'  to  obtain  the 
mir.iir.um  i  value  .  rrn/r>  and  maximum  i  value  i  :fX>, .  for  each  of  the  blocked 

nixies 


3.1  7  ANGLES 

In  'sH  K  1)0  the  hexagonal  mesh  is  skewed  by  60  degrees  to  become  a  carte 
sian  grid  svstem  'scanning  to  define  the  object  boundarx  takes  place  over  a 
given  number  of  angles  in  the  x  v  plane,  which  when  skewed,  become  lines 
•  huh  are  no  ...nger  uniformly  spaced  angularly  This  routine  computes  the 
m<  rements  m  x  and  y  to  be  taken  along  each  of  these  lines  The  increments 
*;e  selected  b»-  'Idu ns ke wed  mesh 


Ax,  =  u  9*1  i  v  3  sin  o,  c.*^,  ;  i  10> 

Ay,  -  0  9**  cos  9,  (111 

W  here  j  n  Ui»-  angP  in  the  \  V  plane  measured  from  the  \  axis  and 

n*  rernente.;  Pi  *  4  gk-  V#  is  the  number  of  scan  angles  to  be 


l%e.J  ,n  tpe  »  i  p.ane 


3.2  SHADO  -  Calculate  Geometrical  Ion  Density 

SHADO  calculates  the  presheath  ion  density  using  the  neutral  ion  approx¬ 
imation.  This  approximation  assumes  that  ion  motion  is  perturbed  only 
by  collisions  with  an  absorbing  object.  The  ion  density  accounts  fully  for 
ion  thermal  velocities,  but  neglects  trajectory  bending  by  electric  or  mag¬ 
netic  fields.  The  ion  density  at  a  point  is  obtained  by  integrating  the  ion 
distribution  at  that  point  over  the  range  of  ion  velocities: 
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n,(£)  =  J  /(f,  v)dv 


The  ion  distribution  is  a  function  of  position  x  and  of  the  ion  velocities  v. 
The  ion  distribution  function  is  redefined  as  the  product  of  the  unperturbed 
velocity  distribution  function  and  a  geometrical  factor  as  follows: 


f(x,v)  =  £  /cum(9)/JV^ 


4  V1 
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The  fcum(9)  array  was  computed  in  FCAL  while  the  geometrical  factor 
is  determined  by  scanning  along  lines  of  constant  <t>  in  the  x’,y’  plane  and 
finding  the  points  at  which  an  unblocked  position  becomes  blocked  or  a 
blocked  position  becomes  unblocked.  If  the  given  point  is  behind  the  object, 
a  6  can  be  found  by  interpolating  for  the  object  boundary  between  blocked 
and  unblocked  scans  and  finding  the  zexit  for  this  interpolated  position. 
The  angle  6  =  arctan(  j^-)  where  r  is  the  distance  from  the  interpolated 
boundary  position  to  the  given  point. 

The  structure  diagram  for  SHADO  is  shown  in  Figure  3. 

The  majority  of  coding  in  SHADO  is  devoted  to  determining  the  status 
of  the  reference  point  and  the  scan  points  on  the  hexagonal  projection  mesh. 
The  status  includes  determinations  of  whether  a  particular  point  is  within 
the  hexagonal  mesh  or  off  of  it;  whether  it  is  in  front  of  the  object;  whether 
the  point  is  blocked  in  the  ram  direction,  and  if  blocked,  whether  the  point 
is  inside  the  object.  As  indicated  in  the  structure  diagram,  the  reference 
point  (input)  is  first  checked.  If  the  point  is  inside  the  object,  in  front  of 
it,  or  so  far  off  to  the  side  of  the  object  that  it  would  not  be  possible  to 
intersect  the  object  within  the  99.9%  limit,  then  the  GI  value  is  returned  as 


Figure  3:  Structure  Diagram  for  the  SHADO  subroutine 
zero  without  further  calculation. 

Once  the  reference  point  status  is  determined,  a  loop  is  entered  which 
initiates  scans  along  the  scan  angles  <t>.  calculated  in  ANGLES.  For  each 
4>i ,  a  transition  is  sought  which  marks  the  outline  of  the  object  on  the 
projection  surface.  If  the  reference  point  is  not  blocked  with  respect  to  the 
ram  direction,  then  a  transition  is  sought  where  the  object  first  blocks  a 
scan  point.  If  the  reference  point  is  blocked,  then  a  transition  to  unblocked 
status  is  sought. 

When  a  transition  is  found,  the  angle  6  is  computed  which  is  the  angle 
from  the  reference  point  to  the  object  surface  at  the  transition  point.  If  the 
reference  point  is  blocked  and  behind  the  object,  the  density  contribution 
for  <t>i  is  fcum{6)l N#.  If  the  reference  point  is  unblocked,  then  density 
contribution  is  1.  -  fcum(0 1)  +  /cum(02)  where  9\  is  the  angular  position 
of  the  first  transition  to  blocked  status,  and  02  is  the  transition  back  to 
unblocked  status. 
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e  7:  Quasisphere  Wake  results  at  z/R=5  using  the  SHADO  algorithm 
N*  =  8 


8:  Quasisphere  Wake  results  at  z/R=5  using  the  SHADO  algorithr 
*  =  16 


['he  shuttle  object  and  its  projection  onto  the  skewed  hexagonal 
ram  direction  is  pointing  directly  into  the  cargo  bay  doors. 
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Figure  12:  The  particle  wake  behind  the  shuttle  using  the  SHADO  routine 
with  N^ex  =  1000  and  N+  -  16 
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